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Abstract: An impulsive model of augmentative biological control consisting of a general 
continuous predator-prey model in ordinary differential equations augmented by a discrete 
part describing periodic introductions of predators is considered. It is shown that there exists 
an invariant periodic solution that corresponds to prey eradication and a condition ensuring 
its global asymptotic stability is given. An optimisation problem related to the preemptive 
use of augmentative biological control is then considered. It is assumed that the per time 
unit budget of biological control {i.e. the number of predators to be released) is fixed and 
the best deployment of this budget is sought after in terms of release frequency. The cost 
function to be minimised is the time taken to reduce an unforeseen prey (pest) invasion under 
some harmless level. The analysis shows that the optimisation problem admits a countable 
infinite number of solutions. An argumentation considering the required robustness of the 
optimisation result is then conducted and it is shown that the best deployment is to use as 
frequent (and thus as small) as possible predator introductions. 
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Stabilite globale et optimisation d'un modele 
impulsionnel general de lutte biologique 

Resume : Un modele impulsionnel de lutte biologique est considere; il est compose d'un 
modele proie-predateur general en equations differentielles ordinaires complete par une partie 
discrete decrivant I'introduction periodique de predateurs. II est demontre qu'il existe une 
solution periodique invariante qui correspond a I'eradication des proies et une condition qui 
permet d'assurer sa stabilite globale asymptotique est donnee. Un probleme d'optimisation 
lie a I'utilisation preventive de la lutte biologique est alors considere. II est suppose que 
le budget par unite de temps de la lutte biologique (c'est-a-dire le nombre de predateurs 
Iciches) est fixe et la meilleure repartition de ce budget est recherchee en termes de frequence 
d'introduction. La fonction de cout a minimiser est le temps pris pour reduire une invasion 
imprevue de ravageurs (les proies) jusqu'a un niveau inoffensif. L'analyse montre que le 
probleme d'optimisation considere admet un nombre infini de solutions. La prise en compte 
de la robustesse necessaire de la solution optimale indique que la meilleure strategie consiste 
a utiliser des lachers de predateurs aussi frequents (et done aussi faibles) que possible. 

Mots-cles : Dynamique proie-predateur, modele impulsionnel, stabilite globale, optimi- 
sation 



Global stability and optimisation of a general impulsive biological control model 3 



1 Introduction 

The biological control method intends to reduce harmful organisms populations and rehes 
on the utiHsation of their natural enemies, the biological control agents. It can be used to 
control e.g. vectors in vector-borne diseases Hke malaria or invasive species (may they be 
animals or plants), but it is mostly dedicated to the control of insect pests, especially at 
the agricultural cropping system scale [25]. Three main biological control strategies can be 
identified: 

• conservation biological control: biocontrol agents are already present in the system 
and are favoured through various means (habitat management...). 

• classical biological control: a small population of biocontrol agents is introduced once in 
the system and a long term equilibrium between the pests and themselves is targetted. 

• augmentative biological control: biocontrol agents are repeatedly introduced in order 
to eradicate the pest population. 

In this contribution we focus on this last method, which is fairly well modelled by an 
impulsive system of ordinary differential equations (a continuous system affected by sudden 
state modifications at some discrete time moments). The intrinsic dynamics of the two 
interacting populations (pests and biocontrol agents) are described by the continuous part 
of the system and the repeated (periodic) biocontrol agents introductions by the discrete 
one. Early examples of the use of impulsive models in biology are [H [12] in chemostat 
modelling. Since then, they have received considerable attention in the area of epidemiology 
[H [ini [28], cancer treatment [6l [3 [17], populations dynamics |20[ [32] or integrated pest 
management |19 [ [27 [ [33 l l37j. to cite some. Though there are numerous studies on bio-inspired 
impulsive models, to our knowledge only a few of these is concerned with optimisation: see 
[H [HI 1211 [HU [3l] . Here we focus both on a global convergence result for a rather general 
predator-prey model with periodic impulsive predator introductions and on a problem of 
optimisation on how to best deploy these introductions in the context of crop protection 
through biological control. Our results follow and generalise earlier ones that were obtained 
in a more restricted context [23^. 

This paper is organised as follows. First, the impulsive predator-prey model of augmen- 
tative biological control that will be payed attention to throughout the paper is presented. 
The hypotheses assumed on the biological part (continuous system) are qualitative so that 
they can encompass a large number of classical biological/ecological functions. The existence 
of an invariant periodic prey-eradication solution is proved and a sufficient condition that 
ensures its global asymptotic stability is given. Then, we focus on an optimisation prob- 
lem in the perspective of preemptive/preventative use of augmentative biological control 
which particularly suits the protection of high valued agricultural crops like e.g. orchards, 
vegetables, mushrooms or ornamentals [U [Til [TH [151 [TBI [29] . It is assumed that a fixed bud- 
get ensuring pest eradication is allowed and that only its deployment (a tradeoff between 
frequent-small or rare-large releases) may be modified. The cost function to be minimised 
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is the time taken to reduce an unforeseen prey (pest) invasion under some harmless level 
and the input is the impulsive release frequency. It is first shown that the local optimisa- 
tion problem has a countable infinite number of solutions and, in a second step using an 
argumentation based on the needed robustness of our result, it is proved that the locally 
optimum solution is to choose as frequent (and thus as small) as possible releases. It is also 
shown that this strategy gives a sub-optimal solution of the global optimisation problem. 
Moreover, the larger the frequency, the lower the cost function, so that even if one cannot 
achieve a very large frequency, it is always interesting to choose the largest possible one. A 
Monte Carlo hke numerical simulation ihustrates our analytical optimisation result and the 
article is concluded by some recommendations to biological control practitioners that follow 
from our analysis. 

2 Model Description and Analysis 

2.1 A general impulsive biological control model 

In this paper, we assume that the tri-trophic system (crop - pests - biocontrol agents) may 
be well represented by its bi-trophic simplification representing the prey (pest) predator 
(biocontrol agent) interactions. The underlying assumption is that the prey population 
remains at moderate levels so that the crop is not limiting and thus do not need to be taken 
into account in the model. We consider the simple case in which no inter specific interactions 
within the predator population affects the predator prey interaction, though some in the prey 
population may occur (see Remark 1 below) . According to classical predator-prey modelling, 
we get the fohowing two-dimensional model: 

( x = f{x) -g{x)y 
\ y = h{x)y -my 

X denoting the preys and y the predators, /(x) denotes the preys' growth velocity, g{x), h{x) 
and m denotes the predators' functional response, numerical response and mortality rate, 
respectively. Since biological processes are always difficult to model, we only assume general 
qualitative hypotheses on the functions /(.), g{.) and h{.): 

Hypotheses 1 (HI): Let /(.), g{.) and h{.) he locally Lipschitz functions on 1R+ such 
that: 

(i) /(0) = 

(ii) g{0) = 0, g'{0) > and Va; > 0, g{x) > 

fix) 

(Hi) the function g^^^ upper bounded for all positive x 
(iv) h{0) = and Va; > 0, h{x) > 
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Remark 1: A large part of the predator-prey functions encountered in the literature fit 
these hypotheses: (Hl-i) only indicates that no spontaneous generation of pests (prey) is 
possible; /(.) might be constant, Hnear, logistic or model an Allee effect eic...(Hl-ii) means 
that there is no prey consumption as the prey level is zero and that if some prey are present, 
the predator are able to find and consume them. Moreover, g{x) is supposed to be increasing 
at the origin. As a consequence g{x) may be a Holling I or II function. g{x) might even 
be non-monotonic (in case of e.g. prey group defence abilities) Hke in Holling IV, but in 
this case from (Hl-iii) f{x) must be of type IV too or become negative for high x values 
(Hke for logistic growth). Notice however that Holling III Hke functional response {i.e. those 
with null derivative at x = 0) fall beyond the scope of this study. No hypotheses are made 
on the numerical response h{x) except (Hl-iv) which is quite natural. Indeed, as it will be 
clearer in the following, our argumentation is mainly based on the fact that g{x) is positive 
for positive x. 

We now model the augmentative biological control procedure as the T-periodic releases 
of biological control agents (with T > 0). Releases are, by their very nature, discrete 
phenomena: every time nT (n G N) some predators are instantly added into the system. 
Let us suppose that there is a fixed rate of predators ^ {i.e. number of predators per unit 
time) to be released in the system. Notice that such a quantity ^ is a direct measure of the 
costs of the biological control method over some time period. Hence in the following it will 
be referred to as well as the release rate or the biological control budget. 

This results in the following: at each time period T, predators are added to the 
predator population y, yielding the discrete process: 



Combining |(T|) and ([2]) yields the following impulsive biological control model, similar 
to, but more general than the ones presented e.g. in |18[ 119): 



Remark 2: It is quite easy to check from (HI) that the proposed model is, as required 
for biological models in general, a non-negative system (see e.g. |22]), i.e. the non-negative 
orthant of the state space is positively invariant by system ^ . Then every non-negative 
initial conditions xq and ?/o at initial time t^ (that only make sense with respect to the 
considered problem) produce through model ^ non-negative trajectories x{t) and y{t) for 
all positive time. 

2.2 Global convergence 

In this section we show that provided the predator release rate ^ is greater than some value, 
augmentative biological control is able to drive any pest population to zero. This result is 
summarised into the following Theorem. 



Vn e N, y{nT+) = y{nT) + 



(2) 




(3) 
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Theorem 1 Under Hypotheses 1, model ^ possesses a "pest free" T-periodic solution: 

/ ;/Tp""(* modT)\ 

{x,{t),y,{t)) = fo / ^ 1 (4) 

which is locally asymptotically stable if and only if: 



and globally asymptotically stable if: 



^ e A rnf{x) 
fi > b — sup — -— (6) 

x>o g{x) 

Remark 3: Note that from (Hl-iii) S as defined in ^ does exist. 

Proof: We first focus on the "pest free" set: {{x, y) e a; = 0}, which is clearly invariant 
by system ^ through (Hl-i) and (Hl-ii). Within this set, according to (Hl-iv), system ^ 
becomes: 

' i = 

y = —my (7) 
Vn e N, y{nT+) = y{nT) + 

that yields: y{{n+l)T^) — y{nT'^)e^"^'^+fiT. It is clear that the sequence (?/(nr+))„gN has 
a single and globally stable equilibrium y* — — e"™"^). Then, system ([7]) trajectories 

globally converges towards the T-periodic solution: 

Vie(nr,(n+1)T], {x{t),y{t))^ O/ 



which is indeed the same as so that the existence of the T-periodic pest free solution 
for system (JS]) has been proved. 

We now concentrate on the stability of the pest free solution ^ for system ([3]) , that is 
to say we do not restrict ourselves to the pest free set. We first change system ^ variables 
to consider the deviations from the pest free solution that are denoted (i, y) so that: 

that yields: 

i = fix) ~ 9ix){y + ypit)) 
y = h{x){yp{t) + y) - my 

We first investigate local asymptotic stability (LAS) of the periodic solution {xp(t), yp{t)). 
Assuming (x, y) are small, we get from ^ at first order in x and y: 

i={f{0)-g'{0)yp{t))x 

y = h'{0)yp{t)x -my ^ ' 
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which is a linear system (in x and y) with T-periodic coefficients. The LAS of {xp{t), Vpit)) 
for the original system ([3]) is equivalent to the LAS of (0,0) for system ([9]). The classical 
approach (Floquet Theory) is to compute the discrete dynamical system that maps the 
variables at some time to the variables one period of the periodic coefficients later. The 
origin is then asymptotically stable if and only if the obtained discrete dynamical system is 
asymptotically stable (see e.g. [l] for the theory and |19| for an application comparable to 
ours). Here, we get: 



((n + l)T+) 




the matrix being lower triangular, the f term does not influence the system's stability so there 
is no need to compute it. Clearly e~™"^ belongs to (0,1), and 

e/o^(/'(o)-9'(o)yp(r))dr -g positive. Then, (0,0) is LAS for ^ {i.e. ixp,yp) is LAS for (0])) iff 
this latter coefficient is lower than one i.e. : 



(/'(O) - g'{0)yp{T))dT < ^ f'{0)T < / g'{0)-^-—^ dr 

^"^ (10) 
to/'(0) ^ ' 

so that the LAS of {xp, yp) part of the Theorem is shown. 

We now focus on the global asymptotic stability (GAS) of the pest free solution. From 
now on, we assume that system ^ is initiated at {xo,yo) at time tg > 0, i.e. system ^ is 
initiated at (xQ^yo) = {xo,yo — yp{to)) at time tg > 0. 
Let us consider the following function: 



r 1 

G(x) = m ——^ ds 



for some k > 0. Since from (Hl-ii) g{x) is positive for positive x, G{.) is an increasing 
function from x = 0, where it goes to — oo since g{.) is locally Lipschitz on ]R+. In the 
following we investigate G{x{t)) behaviour as time t goes to infinity. We have: 

/•£(*) 1 

G{x{t)) - G(io) = m / ^-r ds 

Jio 9[s) 



to 5(i(T)) 
mf{i{T)) 



9{S:{t)) 

Considering system ^ together with (Hl-iii) and Remark 2, we get: 

y > —my 



m{y{T) + yp{T)) ] dr (11) 
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so that, through standard arguments: 

y{t) > min(0,yo)e""(*-*«) 
Using Ijlip and the definition of S" in ([6]) , we then have for all t > to: 

G{x{t)) - G{xo) < J (s'-mzjp(T) -mmin(0,yo)e-™(^-*°)) 



(12) 



{S - myp{T)) dr 



t 




to 




f 







{S - myp{T))dT 



+ (5 - mypir)) dr + min(0, yo) fe^^^*"*") - l) 



(13) 



since yp{t) is T-periodic. It is clear that the first, third and fourth term of the right hand 
side of (fTSl are upper bounded. Now, suppose that ^ holds, then in a similar way than in 

m- 

i-T 

{S — myp{T)) dr < 

so that, since [^J goes to infinity as t does, the right hand side of l|13p goes to — c». Then, 
condition ^ impHes: 

lim G{x{t)) = -oo ^ lim x{t) ^ 

t — *+oo t — > + co 

One can now prove that y converges to zero as well. Indeed, from (fT2|) : 

yo<0 ^ y(t) > yoe-"(*-*«) 

so that, either y{t) converges to zero in infinite time from below, or it reaches the positively 
invariant region where y > in finite time. Up to initial time to translation, we then have to 
consider the positive yo only. Since x converges to zero and h{0) — 0, it is clear that there 
exists a time tf such that: 



Vt > tf, h{i{t)) < 



m , 



yt >tf, y< h{i)yp{t) - —y 



Since h[x)yp{t) goes to zero as t goes to infinity, so does y. We have shown that if condition 
lH]) holds, (0,0) is globally attractive for ([8]), i.e. {xp{t),yp{t)) is globally attractive for 
It is then easily checked through the Hospital rule, (Hl-i) and (Hl-ii) that: 



lim 
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Then, condition Q implies that the necessary and sufficient condition |(5]) for LAS of {xp, yp) 
holds true, which concludes the proof. □ 

Theorem[T]ensures that when the condition ([6]) is verified, the extinction of the prey (pest) 
population is GAS on the positive orthant of the state space. If the local condition |[5]) were 
not satisfied, the pest free solution would then no more be stable. Using bifurcation theory 
for impulsive systems proposed by [T^, one then should be able to perform a bifurcation 
analysis, similar to what is done e.g. in [l9l [26] to show the rich dynamics that system 
^ would produce. However, as stated in the introduction, with respect to our application 
we concentrate here on the pest eradication problem only. In the case where condition ([5]) 
holds true, but condition ([6]) does not, the pest free solution is still locally stable but, since 
our global stability condition ^ is only sufficient, we cannot rule out the possibility that 
it is also globally stable. Such a release rate has the advantage of being smaller than the 
one required for GAS so that it would be cheaper for the grower to choose a one as such: 
it allows good control of limited pests invasions. The price to pay is however the risk of 
considerable crop damages if a large pest invasion occurs. 

It is to be noted that the stability conditions ^ and ^ on the release rate /i are both 
independent of the release period T. Therefore, once a release rate /i fulfilling a stability 
condition is chosen, pest eradication could be achieved independently of the choice of the 
release period T . Indeed, both infrequent and large releases (T and /iT large) or frequent 
and small releases (T and fiT small) would ultimately drive the pest population to zero. 
Though there is no difference in the long term, different values of the release period T might 
however result in different outcomes for the transient dynamics of the system. Hence in the 
following, we will take advantage of the adjustable parameter T to address the practically 
important question: how to "best" (in a sense that will be detailed in the sequel) deploy a 
given biological control budget jj, to eradicate pest outbreaks ? 

3 Optimal Choice of the Release Period in the Preemp- 
tive Case 

3.1 Statement of the optimisation problem 

Throughout the following, we assume that the condition for stability of the pest free solution 
holds true. We focus on preemptive/preventative use of biological control agents releases: 
we suppose that biological control agents are released in anticipation of pests outbreaks, so 
that natural enemies would be able to fight the pest right at the time of their invasion. This 
departs from most classical biological control procedures which are of a more "feedback" type: 
biocontrol agents are released once the pests have been detected only. Such a preventative 
approach has however been payed attention to since it appears as achieving more acceptable 
pest control, especially for high valued crops that are very sensitive to the slightest pest 
outbreak (see [HI |29l SD] for theoretical/simulatory studies and [HI [TH [TH [21 [H] for real 
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life experiments and a biological perspective). Hence in the following, we suppose that 
preventative releases are being performed and that the system is in the invariant set {y > 0} 
as a pest population xq invades the crop at some unforeseen moment io- 

Our goal is to provide recommendations on how to deploy {i.e. how to choose T) such 
a preemptive biological control strategy to minimise crop damages due to a pest outbreak. 
To evaluate crop damages induced by the pests, we use the the concept of "Economic Injury 
Level" (EIL) that has been introduced from the early bases of theoretical biological control 
[30] , EIL (denoted x in the sequel) is defined as the "lowest (positive) pest population level 
that will cause economic losses on the crop", x is then a constant parameter of our problem. 
For a given pest outbreak xq > x, we consider that the damage cost J is an increasing 
function of the time spent by the pest population above the EIL x, i.e. : 

r*f 

J = / ^{T)dT 
J to 

with tf defined such that x{tf) = x and 7(.) > 0. It is clear that minimising J is equivalent 
to minimising 11 = {tf — to) so that we will only concentrate on this latter in the following. 

We are left with two more problems. The first one is that, for given xq and T, the 
damage cost J (or 11) strongly depends on the pest outbreak instant tg, i.e. for a given 
T, the same invading pest level xq yields different damage cost depending on the value of 
to- Without loss of generality, from now on we will consider that S i^^T). We seek to 
provide the safest biological control procedure so that we will concentrate on the damage 
cost for the worst tg, ensuring any other to would result in smaller damage costs. Hence in 
the following, we seek the T that minimises the damage cost for its worst to i. e. we seek to 
minimise max H. 

The remaining problem is that our analysis cannot be achieved for the general model ^ . 
It is however possible to carry it out for a local approximation of the x part of the model 
about a: = as well as for an upper-bound of x. Indeed, consider the local approximation 
of X near {xp, yp). From ([9]), and considering the change of variable zi — In (|) (which 
is obviously increasing in x), we get: 

zi = Si - niyp{t) 

where Si = < S. 

Now we consider the general x equation of model ^ under the hypothesis that the 
system is in the positively invariant set {y > 0}. We get for all t >to: 

i < f{x) - g{x)yp{t) 

Let us make the change of variable Z2 = m 'glxj'^^ (which is increasing in x), we get: 

Z2 < S - myp{t) 
whose right hand side is similar to the ii equation. 
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In the sequel we will then focus on the system: 

i = cr - myp{t) (14) 

with z standing for zi or Z2 and a for Si or the original S. Through equation we will 
thus study in one step the exact local approximation of our system as well as an upper 
bound of it, yielding a locally optimal result as well as a globally sup-optimal one. 



3.2 Preliminary results 

To state our main optimisation result, we first need some preliminary computations. Since 
we assume that the stability condition of the pest free solution holds true, we have: 



Hypotheses 2 (H2): The biological control budget is chosen such that: /i > cr 

We first show that provided the release period T is not too large, the pest population is 
decreasing for t >to. 

Proposition 1 There exists a release period T such that if T < T, the pest population is 
decreasing for t >to. 



Proof: Consider equation Provided mint Upit) > —, z is decreasing for t > t^; then 

so does the pest population x. Notice that: 

is a decreasing function of T with: 

lim mmypit) — — > — 
T-»o+ * mm 

since the stability condition holds. 

Then Proposition [T] holds true with T solution of: 

" (15) 



□ 



Proposition [T] has practical interest regarding growers concerns: indeed it is not interest- 
ing to invest in preventative releases of biological control agents and allow an invading pest 
population to proliferate after the moment of invasion. It is much more desirable to choose 
a procedure that will always make the pest population decreasing after the invasion. Thus 
in the sequel we assume that: 
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Hypotheses 3 (H3): The release period is chosen such that T g (0,T). 

We now seek the moment of invasion to that maximises the damage costs J {i.e. that 
maximises the damage time 11) for a given invading pest population xq > x (corresponding 
to zq > 0) and a given release period T e (0, T). Notice that with the z formalism, we have 
li — t f — with z{tQ) = zq and z{tf) — 0. We have: 

Lemma 1 [[23]] Suppose zq > and T G (0, T) are fixed, then one of the following holds: 

(i) 3k e N, (^ - (j)kT ^ zo then H = fcT 

(ii) 3fc e N, maxto n(to) = n(t*) = (fc + 1)T - t* («.e. z((fc + 1)T) = 0) 

We will first show case (i) and prove that, otherwise, either case (ii) holds or 11 is 
maximum at t*^ — 0, this latter case being impossible. 



Proof: Suppose that n(io) is maximum for tg G (0, T) and such that [t^ + n(tQ)) S 
(fcT, [k + 1)T) for some integer k. That is to say that we suppose that n(to) is maximum 
for both io and tf — (n(to) +to) being strictly within two moments that are multiples of T. 

Now pick <tl< t^^ in (0,r) such that (to + n(to)) G [fcT, [k + 1)T) for all to 6 
[f^^t^']. Within this set, define ^(to) as: 

n(to) = fcT - to + 0(to) (16) 

Integrating (fT4| between to and to + n(to), we get: 

rta+n(to) 

z(to + n(to)) = Zo + (jIl{to) - m yp{T)dT 



Zo + crn(to) 



1 



f e-™^dT + (fc-l) f e-'^^dr^ 

Jto Jo ^ 

^n(to) - YZ^i^ (^"""'' - e"™'^'"^ + - ^""^)) (17) 



which is, from the definition of II, equal to 0. 
From i|16p . we have: 

dH , . dO , . ^ , . 

Using this and differentiating (fTTjl with respect to to yields: 

^(io) = -r (19) 

dto e~™»(*o) - cr 
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To have a maximum of 11 at ip, it is required that: 

so that from ifTO]) . we must have 0{tQ) — tg. Then, from lfT7)l and ifTG]). the integer A: must 
be such that: 

(/i — cr)fcr = Zq 

and 11 = kT does not depend on to, which show case (i). Otherwise n(to) has no extremum 
at t^ verifying t^ £ (0,T) and {t^ + n{t^)) e {kT, (fc + l)r). 

Two cases remain to be studied: either 11 is maximum for to = 0, or t^ is such that 
(^Q + 11(^0)) = (fc + l)T. The remaining cases t'^=T and (tQ + n(<g)) — kT might be studied 
by k reparametrisation. 

We concentrate first on the case = 0. Then from (fTSl and ifTOl) the right derivative of 
n at to = is: 

/ m^T \ 

— (to = 0+) = ^ ^ ^ 1 (20) 

From Proposition [ll we have VT G (0,r): 

_mfjJ_ J, 
a < 7- 

(1 - e-™^) 

Since by definition 6'(to) < T both the numerator and denominator of the fraction in the 
right hand side of (|20l) are positive, the former being larger than the latter, ^(^o = O"*") is 
then positive so that 11 has a (local) minimum at to = 0, ruling out the case tg = 0. 

We now focus on the case where tg is such that tg + n(tQ) = (fc + 1)T. Then from l|16p . 
0(to) = T and the left derivative of 11 at tg is: 



( m/jT ^ 






( in/iT ^ 


1 e-™^ - (J 





In a very same way as in the previous case, since T < T we show that this derivative is 
positive. Similarly one can show that the right derivative of 11 at tjj is negative. This can be 
performed through the reparametrisation of A: as {k + 1) while noticing that this corresponds 

to e(t^) = 0. 

Then for the considered zo and T, 11 is maximum for t^ such that (t^ +n(t5)) = {k + l)T, 
which completes the proof of case (ii). □ 



Lemma[T]is quite natural. Indeed, from the linearity and periodicity of ifH]). the decrease 
that takes place between times T and kT is independent of to. The worst case should then 
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contain the end of the first time interval where the predators are scarce, rather than the 
beginning of the last interval where they are abundant. 

In the following we look for the release period T G (0, T) that minimises the worst case 
damage time: maxt^ 11. 

3.3 Optimal Choice of the Release Period 

We now investigate which T G (0,r) minimises the worst case damage time maxj^ IT {i.e. 
that minimises the worst case damages maxt^ J). We get the following result: 

Theorem 2 Suppose zq > is fixed. Let: 

Ti = (21) 

fi — a 

T 

Then there exists uq G N* such that min max n = Ti is reached at T = r„ = — for all 

Te(o,f) to 72 

integer n > tiq. 

Proof: Consider Ti that is positive since /i is assumed to be greater than a. It is then 
clear that there exists an integer no such that for all integer n > nQ,Tn <T. Now consider 
a T„ with n > uq; using (|2T1) . we have: 



in - a)nTn = (^i - <j)Ti = zq 

Such a Tn corresponds to case (i) in Lemma [1] thus we have 11 = Ti which does not depend 
on to- 

Now we show that 11 = Ti is the minimum (with respect to T G (0, Tj) of maxt^ U. Pick 
aT ^ Tn for all n > uq. For such a T, case (i) of Lemma [T] is ruled out and according to 
case (ii), we have: 

3k G N, maxn = U{t*) = {k + 1)T ~ t* 

to 

and z{{k + 1)T) = 0. Using ^ and 6'(i5) = T, we have: 



.1 - e 



zo + cT {{k + 1)T - t*) - (fc + l)^iT + ^Tj—-^ = (22) 



so that: 



which yields: 



(fc + 1)T = [zo-c7t*+fiT 



,1 - e- 



maxn = {k + l)T-t* 

to 



= ^^ + ;7^^T3^^-*S) (23) 
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Notice that the function (^^—^ — ^ is decreasing of t, so that, since < T, we have: 

-mt' \ 



1 - e- 

Since from (H2), /i > ct, we have shown that max^g 11 is greater than Ti for all period T 
such that \fn > no, T ^ Tn, which completes the proof. □ 



Theorem [2] shows that for a given zq we have a countable infinite number of release 
periods, the T„, that solves our optimisation problem. 

However, we have two difficulties with Theorem [21 On the one hand, through equation 
lfT4|) we studied both the local optimisation problem and the global suboptimal one. To 
solve these two, we need to choose a release period T = T„, which depends on a. In the 
local optimisation case, a equals Si = "I^q^-* while in the global suboptimal one, a equals 

S — sup2,>Q . There is no reason why would be rational, so that Theorem [2] states 
that it is not possible to choose a release period T that solves both the local and global 
optimisation problem (what remains however the objective). 

On the other hand, it is to be noticed that Ti, which determines the optimal values of 
the release period, depends on some "biological" values, namely the invading pest population 
level xq {i.e. zq) and the combination of parameters a. This is quite a big drawback since 
biological parameters are usually badly known and the initial invading pest level a:o is not 
known at all. Clearly such model-based optimisation approaches should take into account 
these uncertainties [2Tj . 

Hence, in the following, we will show how T should be chosen in order to ensure that the 
worst case damage time maxt^ U remains close to its minimal value, independently of the 
value of Xq and a. In doing so, we will provide an almost optimal choice of T that solves 
both the local and global optimisation problem and that is robust to the value of xq and to 
parameters uncertainty. 



3.4 Robustness Analysis 

In the preceding sections, we wanted our choice of T to minimise the difference between 
the time that it takes to reach x from xo for the worst possible choice of the initial time 
to. In addition to that, we now also want this to be robust with respect to the uncertainty 
on the initial condition xo and to the model parameters that we wih denote p in the sequel 
{p is simply the pair (cr, m). Indeed, in the poorly measured conditions of crop culture and 
with the accompanying scarcely known models, it is fundamental for a control method to 
be robust with respect to these uncertainties. We then suppose that the "true" a;o and p 
belongs to some compact sets: 



Hypotheses 4 (H4) : The uncertainties on the initial condition and the model parameters 
are bounded, with zo G [zo, zq] denoted Z and p that belongs to a compact set V . 
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We of course still assume that the invading population is greater than the EIL x {i.e. 
£0 > 0) and that the stability condition holds, i.e. /i > maxpg-p (ct) . Our optimisation 
procedure will be robust provided our choice of T minimises maxt^ 11 for the worst zq and p 
choice within the sets defined in Hypotheses (H4). 

Notice that the minimum (with respect to T) of max^^ 11 is equal to Ti that does depend 
on both zq and p. Thus in our robustness analysis it makes more sense to consider the 
deviation of maxt^ 11 from Ti, rather than the absolute value of maxt„ 11. Hence, we first 
focus on the evolution of: 



max 



(^maxU{to,zo,p) ] -Ti{zo,p) 



with respect to T. We have: 
Theorem 3 The function: 



max 



(ma.xIl{to,zo,p)j ~~Ti{zo,p) 



is increasing for T smaller than some positive constant T^. Moreover: 



lim max 

T^0+ {zo,p)eZxV 



ma.xTl{to,ZQ,p) ] -~Ti{zq,p) 



(24) 



(25) 



Proof: From equation l[23|) . we have: 



max 



(^niaxn(to, zo,p) ) -Ti(zo,p) 



max 

{zo.,p) ezxv 



/i — (T 



1 



T - t^{T, zq,p) 



max 



We first concentrate on: 



max 

zo£Z 



^ — a zo e z 

I _ p-rntl(T,za,p) 



-mt'g{T,zo:P) 



^ max { - — ^ T - tQ{T, zo,p) 



T -tl{T, za,p) 



(26) 



for fixed T and p. Let us introduce: 



Ti{zo,zo,p) = min T{p) 



zq- Zq 
2(m - ^) 



with T{p) defined in equation l|15p . In the sequel, we consider that T belongs to (O, Ti(zo, £o,p)) . 
On the one side, this ensures through Proposition[T]that, independently of the actual param- 
eters, the pest population {i.e. z in our new variables) always decreases after the invasion 
time to- 
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On the other side, we have: 

Then, there exists n G N and zqi, ^02 G [£o,^] such that: 

n and r- = n + 1 



T(/z - a) Tip - a) 

From Theorem [21 T is thus an optimal period for both the initial conditions 2:01 and zq2- 

We first show by contradiction that for all zq e (2:01, Z02), case (ii) of Lemma [T] holds for 
k = n. Indeed suppose that for a given zq G (zqi, 2:02), we have k < {n — 1), then: 

maxU = {k + 1)T - t*„ < (n - 1 + 1)T - 

*0 

< Ti(zoi,p)-tl 

< T,{zo,p)~t* (27) 

since from Ti is an increasing function of zq. |(27|) is not possible since from Theorem 
[21 maxll > Ti{zo,p) and t^ > 0. Therefore k cannot be smaller or equal to (n — 1) for 

to 

ZO G {zqi,Z02). 

Now suppose that k > {n+ 1) for some zq € (-Zoi, •202)- From ((22l) . we have: 

1 _ e^™*i' 

Zo = (/i — o')fcr + (/i — (t)T + atQ — piT 



1 - e-™J 

1 -e- 



> (^-a)(n + l)r-(7(T-t*)+^T(^l 



(28) 



The bracketted expression in l(28ll is positive since: 

^m{T-tl) _ I 



> m 
and from (fTSj) : 

VTe (0,f), mr> -(e™'^-l) 
M 

Hence, l(28l) contradicts the hypothesis that zq < zq2 so that k cannot be greater or equal 
to (n + 1) for Zq G (201,202). Thus, we have shown that for all zq S (2^01, -2^02), case (ii) of 
Lemma [H holds for k — n. 

Since T is an optimal period for the initial condition 201, we should have by continuity: 

lim^ n(to(T, 2o,p)) = (n + l)T- lim^ io(T, 2o,p) =Ti{zoi,p) 
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Then: 



lim t^{T,zo,p) = {n + l)T -Ti{zoi,p) 



= {n + l)T- 

A similar argumentation on zo2 yields: 

lim io(T, zo,p) = 



zqi 
H — a 



From (|22l) . t5(r, zo^p) is a continuous function of zq. Then, as zq covers [201,202], to{T, zo,p) 
reaches its whole interval of definition [0, T] so that we have: 



max 

^0 [2^01 j^02] 



I _ Q-™.tl{T.,za,p) 

^ _ ^_„,r T-i*(T,zo,p) 



max 



1 _ e~'"*o 



-mT ^ *0 



T ^ - tl{T,zo,p) 



max 

zo G 2 



1 - e- 



Differentiating 
for: 



T — to ) ■^ith respect to , we show that it reaches its maximum 



= — In 



mT 



1 _g-mT 

il > since mT > 1 - g-™^. In addition, t* < T since e™^ > j^f^. Then: 

' 1 — e^™'c''^-^'^°'^^ 



(29) 



max 

zaeZ 



1 -e- 



1 - e^™**! 



1 

m 



mT 



l_g-mT 

H{T,p) is an increasing function of T on (0, r;(£o7'zo,p)) since: 



- 1 - In 



mT 



1 - e~ 



dH{T,p) e"™^ (e™^ - 1 - mT) 



dT 



(1 - e-^Ty 



1 



1 - e 



-rnT 



mT 



> 



Moreover, through the Hospital rule, we have for all p eV: 

lim H(T,p) = 

0+ 
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Let us introduce: 



Tl = minT;(zo, zo,p) > 
pev — 



Through the above analysis, it has been shown that for T £ {0,Tl) and p & V, H{T,p) is 
an increasing function of T with Iim7-^Q+ H(T,p) = 0. 

We now come back to the function defined in ll24l). We have: 



max 



\^8ixIl{to,zo,p) ] -Ti(zo,p) 



HiT,p) 



Then, from H{T,p) properties, for T smaller than Tl the worst deviation of maxf^ IT from 
Tl (according to {zo,p) € Z x V) is an increasing function of T. Moreover: 



lim max 

T^o+ {zo,p)ezxr 



^niaxn(to, 2o,p) ) -Ti{za,p) 



= 



which concludes the proof. 



□ 



One can deduce the following Corollary from Theorem [3l 

Corollary 1 Let T = minT(p). Then: 

pev 



inf_ max 

Te{o,T){zo,p)eZxv 



max Il{to,zo,p) -Ti{zo,p) 
to e [o,T] 



is achieved for T = 0+. 
Proof: From Theorem [3l 



max Il{to,zo,p) ] -Ti{zq,p) 
to e [o,T] 



inf max 

Teio,TL)izo,p)ezxv 



and is achieved for T = 0+ since the argument of the infimum in l|3T|) is increasin; 
tends to as r does. 

To complete the proof of Corollary [T] we have to study: 



(30) 



(31) 
: of T and 



max 

{zo,p)eZxv 



max Il{tQ,ZQ,p) ] ~Ti{zo,p) 

to e [o,T] 



for T e [Tl,T). Notice that from Tl definition, we have: 

Tr = min ( min ( T(v), — =- I I = min ( T 

(zo,p)eZxv\ \ ^"^"2(11- a) 



' (zo,p)eZxv \2(^ - (t) 
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so that the set [T^.T^ is non-empty if and only if: T > min(^g pjg^rxp (^ 2{fj.-^) ) • Suppose 
that this holds and that T £ [Tl,T), then we have: 

' I _ g-mto(T,zo,p) 



max 

zo£Z 



1 



-mT 



t*oiT,zo,p) >0 



since otherwise it would be required that Vzq G Z there would exist n e N such that 
T = which is obviously not possible. Then, VT e [Tl,^): 



max 

{zo,p)eZxV 



max n(io,2:o,p) -^1(20,^) 
to e [o,T] 



max 



A* 



/X — cr zoe2 

this latter being positive, which concludes the proof. 



-mt^{T,zo,p) 



a 



With the help of Theorem [3l we know that provided T is not too large, the smaller 
T is chosen, the less will be the worst case (according to parameters and initial condition 
uncertainty) deviation of maxt^ 11 from its minimal value Ti . It has moreover been shown 
that as T tends to 0, the worst case deviation of maxt^ 11 from its minimal value goes to 
as well. 

We conclude from this analysis that the robust choice of T, as defined by the criterion 
l(30|l consists in choosing T — 0. However, it is obviously not possible to achieve such a T 
in practical applications: it would consist in continuously releasing predators into the crop 
culture. As an alternative, and since the argument of the infimum in ((3Ql) is an increasing 
function of T for T G (0, T^), it is then advised to take T as small as possible to achieve the 
best possible robustness. 

A last remark should be made about the need for robustness. In the preceding analysis, 
we have concentrated on a single model in the form: 

z = S - yp{t) 

We should remember that this model covers two cases: a local linearisation of the system 
(which yields a certain value of S) and a system that globally upper-bounds the actual 
nonhnear system (which yields another value of S). At the very least, the choice of T 
that we made should work well in both those cases. We can then conclude that, even in 
the unrealistic case where the initial conditions and the parameters are very well-known, 
one should make a choice of T that is sufficiently robust, so that both global and local 
approximations are covered, i.e. one should take T as small as possible. 



4 Numerical Simulation 

To illustrate the robustness result stated in Theorem [3] we performed a Monte Carlo fike 
numerical simulation: an impulsive predator-prey model of the form ([3]) was considered with 
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constant a priori chosen parameters for the continuous part of the model and ^ verifying 
condition A triplet of the remaining parameters was randomly chosen: first a release 
period T G (0, Tl), the pest invasion moment to g [0, T) and the pest invasion level xq > x. 
The value of Tl is easy to compute from (flSl) since in this case the parameter set V is 
reduced to a singleton and that we have chosen a large Z interval. The time to reach x, H, 
was then computed as well as its deviation from Ti . The process has been repeated until a 
set of 2.10^ different simulations was obtained and the results, the deviation of 11 from Ti, 
were plotted against the time period T. 

The results of this numerical work is presented on Figure[ll each of the dot corresponding 
to one of the 2.10^ simulations. What actually corresponds to Theorem [3] is the upper 
envelope of the data set only (since it is the maximum with respect to parameters of the 
worst case, according to io, of 11). 




Figure 1: Graphical illustration of Theorem [3l Each dot corresponds to the deviation of 11 
from the minimal worst case Ti for one of the 2.10'^ randomly chosen triplet T G (0, Tl), <o S 
(0, T) and xq > x. 

As it can be noticed the smaller the choice of T, the smaller the positive deviation of 11 
(thus of the crop damages) from its optimal value Ti . The price to pay is however that such 
a choice reduces also the negative deviation of 11 from Ti . Indeed an analysis similar to the 
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one performed in Section [3l but focusing on: 



nun 




minn(to, zo,p) 



) 



Ti{zq,p) 



would show that this is a decreasing function of T with limit as T tends to 0"*". Hence, if 
one feels lucky, a choice of T small may not be the best one since one can obtain, fortunately, 
smaller 11 (thus smaller crop damages) if the actual is not the worst one. 

5 Conclusion 

In this contribution we studied a rather general predator-prey model with periodic impul- 
sive additions of predators (releases) that is relevant for biological control modelling. The 
impulsive releases of the biological control agents were modelled on the basis of a fixed bud- 
get to be spent per time period (release rate). As such it allowed to compare the different 
releases period as the same overall amount of predators was used. In a first step, it was 
shown that, provided the release rate of biocontrol agents was high enough, the biological 
control program allows to eradicate the prey (pest) population, whatever the value of the 
release period. Then the efficiency of different release periods (for a given release rate) were 
compared on the assumption that the biological control program is conducted on a preven- 
tative basis. It was shown that the time needed to eradicate an unforeseen pest invasion 
occurring at a worst moment had a minimum that could be achieved through the choice 
of a countable infinite number of release period. However, these release period values were 
strongly dependent on the model parameters and a robustness analysis of this optimisation 
problem showed that the smaller the release period, the less the time needed to eradicate 
the pest invasion, for the worst invasion moment and accounting for parameter uncertainty. 
Hence, as an advice to biological control practitioners, our study indicates that it is less risky 
to release frequent small amounts of biocontrol agents rather than massive infrequent ones. 
This result (frequent releases of small doses yield better control) is similar to the ones of 
|35 [ [36 l [TT| that were all obtained from simulatory models. Such a biocontrol strategy is also 
advocated by [2| from real life experiments and by [29] on the basis of a simulatory model. 
However, to our knowledge, this contribution is the first to recommend such a strategy on 
the basis of a mathematical analysis. 

In a broader context, we believe that the techniques used in this work may also be 
of interest for optimisation problems on impulsive models from other applied fields: for 
instance, optimisation of vaccination strategies in epidemiology |28[ [H [13] or optimisation 
of chemotherapy in cancer treatments ^17, .6. .5j. 
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